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ABSTRACT 

A clear prediction of the Cold Dark Matter model is the existence of cuspy dark matter halo density profiles 
on all mass scales. This is not in agreement with the observed rotation curves of spiral galaxies, challenging 
on small scales the otherwise successful CDM paradigm. In this work we employ high resolution cosmological 
hydro-dynamical simulations to study the effects of dissipative processes on the inner distribution of dark matter 
in MiUcy-Way like objects (M IO'^Mq). Our simulations include supernova feedback, and the effects of the 
radiation pressure of massive stars before they explode as supernovae. The increased stellar feedback results in 
the expansion of the dark matter halo instead of contraction with respect to N-body simulations. Baryons are 
able to erase the dark matter cuspy distribution creating a flat, cored, dark matter density profile in the central 
several kpc of a massive Milky-Way like halo. The profile is well fit by a Burkert profile, with fitting parameters 
consistent with the observations. In addition, we obtain flat rotation curves as well as extended, exponential stellar 
disk profiles. While the stellar disk we obtain is still partially too thick to resemble the MW thin disk, this pilot 
study shows that there is enough energy available in the baryonic component to alter the dark matter distribution 
even in massive disc galaxies, providing a possible solution to the long standing problem of cusps vs. cores. 

Subject headings: cosmology: theory — galaxies: structure — hydrodynamics — methods: numerical 



1. INTRODUCTION 

The theory of cold dark matter (CDM) provides a success- 
ful framework for understanding structure formation in the uni- 
verse (e.g. Spergel et al. 2003, Komatsu et al. 2011). In this 
paradigm, dark matter first collapses into small haloes, which 
merge to form progressively larger haloes. Galaxies are thought 
to form out of gas which cools and collapses to the centres 
of these dark matter haloes (White & Rees 1978). Properties 
of dark matter haloes have been extensively studied via colli- 
sionless A^-body simulations of the growth of primordial fluc- 
tuations into gravitationally bound structures. Numerical sim- 
ulations have facilitated detailed predictions for a wide range 
of properties of dark matter haloes at all epochs (e.g. Maccio 
etal. 2008,Prada et al. 2011). 

Dissipationless cosmological simulations have also raised prob- 
lems for the CDM scenario on small scales, one of which is the 
central slope of the dark matter density profile of virialized ob- 
jects. A^-body simulations predict a central concentration, with 
a logarithmic slope of « -1 (Navarro, Frenk & White 1997, 
Diemand et al. 2005, Springel et al. 2008). Such a 'cuspy' 
matter distribution is not supported by observations of the ro- 
tation curves of spiral galaxies, which have revealed that the 
dark halos encompassing disc galaxies have a constant den- 
sity core (e.g. Salucci & Burkert 2000, Oh et al. 2008, Spano 
et al. 2008, Kuzio de Naray et al. 2009). Donato et al. (2009) 
have recently analyzed the rotation curves of a large sample of 
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galaxies ranging over all Hubble types and with luminosity as 
high as Mb ~ -22. Their analysis shows that the central surface 
density, is nearly constant and independent of galaxy luminos- 
ity. This issue presents a major challenge for the otherwise suc- 
cessful CDM cosmological model. 

By construction, dissipationless simulations do not include 
baryons. While on large scales, the effect of gas and stars can 
be neglected, this is not true on small scales, where baryons 
can be gravitationally dominant. For this purpose, cosmologi- 
cal hydro-dynamical simulations have been extensively used to 
directly address the question of galaxy properties in the CDM 
scenario (e.g. Navarro & Steinmetz 2000, Abadi et al. 2003, 
Brook et al. 2004, Robertson et al. 2004, Okamoto et al. 2005, 
Maccio et al. 2006, Governato et al. 2007, Sanchez-Blazquez 
et al. 2009, Scannapieco et al. 2009, Stinson et al. 2010, Pio- 
ntek & Steinmetz 201 1, Agertz et al. 201 1). 

The response of dark matter to baryonic infall (and star for- 
mation) is still highly debated. During galaxy formation, as 
cosmic gas cools and condenses towards the halo centre and 
forms stars, dark matter particles are pulled inward and increase 
their central density. This process is dubbed 'adiabatic con- 
traction' (e.g. Blumenthal et al. 1986, Gnedin et al. 2004). 
Halo contraction is present in the vast majority of cosmolog- 
ical hydro-dynamical simulations (see Gnedin et al. 2011 and 
references therein), yet these simulated galaxies fail to repro- 
duce observed rotation curves due to their too centrally con- 
centrated stellar and dark matter profiles. In order to repro- 
duce observational data, like the rotation velocity-luminosity 
and size-luminosity relations, models with NO adiabatic con- 
traction (Gnedin et al. 2007, Courteau et al. 2007) or even with 
expansion (Dutton et al. 2007, 201 1) are required (but also see 
Trujillo et al. 2010). 

To resolve these discrepencies, theoretical arguments and sim- 
ulations have proposed baryonic processes that can result in 
producing an expansion of the dark matter halo. Gas bulk mo- 
tions, possibly supernova-induced in regions of high star forma- 
tion activity, and the subsequent energy loss of gas clouds due 
to dynamical friction can transfer energy to the central dark- 
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matter component (Navarro, Eke & Steinmetz 1996, Mo and 
Mao 2004, El-Zant et al. 2001,Ogiya & Mori 201 1). 

Mashchenko et al. (2006, 2008) have pointed out that there 
is another (possibly more relevant) effect, namely the gas bulk 
motion can induce substantial gravitational potential fluctua- 
tions and a subsequent reduction in the central dark matter den- 
sity. Recently Governato et al. (2010) have presented high res- 
olution cosmological simulation of a dwarf galaxy that was 
able to create a cored dark matter profile at z = within the 
CDM scenario, and reproduce several properties of observed 
dwarf galaxies. Resolution of an inhomogeneous interstellar 
medium, and strong supernova driven outflows, which inhibits 
the formation of bulges, resulted in the decrease of the dark- 
matter density to less than half of what it would otherwise be 
within the central kiloparsec of these low mass objects (M,;,„ « 
1O"^M0). Further analysis of this simulation, Pontzen & Gov- 
ernato (2011) showed that the flattening was the result of rel- 
atively small starbursts in the centre of the protogalaxy, which 
contribute over many cycles to a gradual transfer of energy from 
the baryons to the dark matter. This mechanism is closely re- 
lated to matter outflows, but does not require violent, sudden 
mass loss. In a recent paper Brook et al. (2011) showed that 
large fractions of the gas that is expelled from the central re- 
gions of galaxies returns via a large scale galactic fountain to 
form stars at later times: this greatly increases the occurrence 
of outflows from the inner regions for a galaxy of given stellar 
mass, and thus any flattening mechanism that relates to outflows 
will be significantly enhanced. 

Dekel & Silk (1986) showed that supernovae can eject mat- 
ter from halos up to 100 km s"', but it has yet to be seen what 
impact this might have on dark matter profiles, nor how the 
addition of radiation pressure feedback might change things. 
So, while observations show evidence for flattened dark matter 
density profiles up to L* galaxies, the question remains whether 
there is enough energy input from baryons in more massive ob- 
jects in order for these processes to be effective in altering the 
dark matter density profile of spiral galaxies with a dark matter 
mass of the order of 10" - 10'^. 

In this letter we present a high resolution cosmological hy- 
drodynamic simulation of a massive spiral galaxy that includes 
cooling due to hydrogen and heavier "metal" elements (Shen 
et al. 2010), UV background radiation (Haardt & Madau from 
CLOUDY), a simple commonly used star formation prescrip- 
tion, adiabatic supernova feedback (Stinson et al. 2006), along 
with feedback from the early radiation produced by massive 
stars (Stinson et al. in prep, see Brook et al. 201 1). Our results 
show that reasonable baryonic feedback is able to create a den- 
sity core in the dark matter distribution even for massive spiral 
galaxies approaching the mass of our own Milky Way. 

This letter is organized as follows in Section |2] we describe 
our simulations, in Section[3]we present our results on the dark 
matter radial density profile, finally, Section|4]is devoted to con- 
clusions and discussion. 

2. SIMULATIONS 

In this paper we use g5664, a cosmological zoom simula- 
tion drawn from the McMaster Unbiased Galaxy Simulations 
(MUGS). See Stinson et al. 2010 for a complete description of 
the creation of the initial conditions. g5664 has a total mass 
of 7 X 10" Mq, a spin parameter of 0.024, and a last major 
merger at z = 3.4. Inside ri,,> at z = 0, there are 4 x 10^ dark mat- 
ter particles with a mass of 1 . 1 x 10^ Mq, a similar number gas 



particles with mass 2.1 x 10^ M0 and between 3.5 x 10^ and 
10^ stars with mass 5.5 x 10^ M0 depending on the star forma- 
tion and feedback recipes used. Using the physics employed in 
the original MUGS simulations, g5664 formed an exponential 
disk with an exponential bulge with a central surface brightness 
of fjLi =18, had a total face-on magnitude of M,-=-2 1 .7 and g-r 
colour of 0.52. In many ways, g5664 is a similar initial condi- 
tion to the Guedes et al. (201 1) eris simulation, except 8 times 
lower resolution. We will refer to this simulation as the low 
feedback run (LFR). 

The new version of g5664 was also evolved using the smoothed 
particle hydrodynamics (SPH) code GASOLINE (Wadsley et al. 2004). 
For this high feedback run (HFR), three changes were made 
which resulted in a stronger implementation of stellar feedback: 

i) The Ki-oupa et al. (1993) IMF used in MUGS was changed to 
the more commonly used Chabrier (2003) IMF which creates 
more massive stars hence more energy per stellar mass created. 

ii) We use a star formation density threshold of n,/, 9.3 cm"^ and 
an efficiency of = 0.1. iii) The new runs assume the energy 
input from supernovae is 10^^ ergs instead of the 4 x 10^" ergs 
used in MUGS, iv) the new run includes radiation released by 
the massive young stars before they explode as supernovae. All 
of these changes are based on a parameter search that will be 
reported in Stinson et al. (in prep) It will be shown that each 
of these changes is necessary to produce more realistic galax- 
ies as described in ^ For details on the implementation of the 
radiative feedback, see Brook et al. (2011), though the work 
presented here deposits 17.5% of the massive star luminosity 
as thermal energy instead of the 10% reported there. Again, 
this is motivated by the production of a more realistic galaxy. 
We note here that the energy deposition from supernovae and 
radiation pressure couple inefficiently to the gas. Both types 
of feedback are deposited into the high density gas found in 
the disk where the cooling time is shorter than the dynamical 
timescales resolved in the simulation. 

We emphasise that our feedback scheme relies on subgrid 
physics. It is not possible to fully resolve the sites where star 
formation happens and energy is fed back into the ISM. Our 
simple description of the radiation pressure mechanism is a rel- 
atively crude initial representation. Finally, we also ran a dark 
matter only version of g5664 (N-body run), with the same par- 
ticle resolution. 

3. RESULTS 

In this Letter we focus on the effects of feedback on the dark 
matter distribution in our simulated galaxy. A more compre- 
hensive study of the properties of simulated galaxies using such 
a strong feedback implementation will be presented in a forth- 
coming series of papers. In brief, the simulation presented here 
does a much better job reproducing the observed properties of 
galaxies including a flatter rotation curve, an exponential sur- 
face brightness profile and a stellar mass, 6 x lO'' Mq, in bet- 
ter agreement with what halo abundance matching predicts a 
7 X 10" Mq should contain (Moster et al. 2010). This latter 
point is important, as it indicates that the large scale outflows 
inherent in this study may be necessary. 

Figure [T] shows the dark matter density profile in our three 
simulations. The pure dark matter run (N-body, black line) 
is well fit by an NFW profile with concentration parameter of 
c = 6, in agreement with cosmological expectations (Munoz- 
Cuartas et al. 201 1). The profile of the low feedback run (LFR, 
dotted-blue) shows evidence of significant adiabatic contrac- 
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Fig. 1. — Density profile for only the dark matter in the three different 
realizations of our galaxy. The blue line shows the low feedback run (LFR). the 
black line shows the dark matter only run (N-body) and the red line shows the 
higher feedback case. The blue curve shows evidence for adiabatic contraction, 
the black one presents the usual NFW profiles, while the red one shows a 
clearly a cored profile, in agreement with observations. 



Q. 



o 
a 




tion, with DM pulled towards the inner regions by the cen- 
trally concentrated baryons. The inner profile is fit with a sin- 
gle power law (p oc r""), with a = 2. As reported in Stinson 
et al. (2010), this dark matter peak is accompanied by a high 
concentration of baryonic material at the centre of the galaxy, 
represented by the a centrally peaked rotation curve and high 
bulge-to- total ratio. None of these features agree with observa- 
tions, which do not support the adiabatic contraction scenario 
at these mass scales. The lowest curve is our HFR, which uses a 
Chabrier IMF and radiation pressure feedback. The dark matter 
density profile follows the pure dark matter run in to r w 5kpc, 
but then it notably flattens to clearly reveal the presence of a 
core in the inner region. 

The DM density profile of the HFR can be fit with a Burkert 
profile (Burkert 1995): 



This profile, when combined with appropriate baryonic gaseous 
and stellar components, is found to reproduce very well the os- 
erved kinematics of disc systems (e.g. Salucci & Burkert 2000; 
Gentile et al. 2007). The two free parameters (/0();ro) can be 
determined through a x' minimization fitting procedure: In our 
case this led to po/pcr = 1-565 x 10^ and tq = 9.11 kpc. The 
simulated DM profile with its Burkert fit are shown in the up- 
per panel of Figure |2] 

Results of Donato et al. (2009), showing the central surface 
density po, defined as the product of the halo core radius and 
central density (po = ropo) of galactic dark matter haloes, are 
shown (open circles) in the lower panel of Figure |2] Our sim- 
ulated galaxy is over-plotted as a red star Not only can the 
simulation be fit with a Burkert profile, but the cored profile 
of g5664 HFR agrees with observed density profiles. The ma- 
genta squares are lower mass simulations which have similar 
high feedback prescriptions as the HFR, with slight calibra- 
tion changes as these simulations have 8 times better resolution. 
Detailed properties of these simulations will be presented in a 
forthcoming paper (Brook et al. 201 lb in prep). 



Fig. 2. — Upper panel: density profile of the DM component in g5664 
HFR, and fitting Burkert profile with a core size of ro = 9.11 kpc. Lower 
panel: The relation between luminosity and dark matter halo surface density. 
Open symbols represent observational results, while our simulated galaxy is 
represented by the red star. The dashed line is the fit to this relation, suggested 
by Donato et al. (2009). The magenta squares are lower mass simulations 
which have similar high feedback prescriptions as the HFR (Brook et al. 201 lb 
in prep). 

Figure [3] shows the dark matter density profile for the hy- 
drodynamical simulations (low and high feedback) at z = 4.8 
and z = 1 (upper and lower panel respectively). The two runs 
show markedly different behaviours: High feedback results in 
low star formation rates in low mass progenitors, as it prevents 
significant gas cooling to the very central regions of the dark 
matter halos. The dark matter profile remains unperturbed from 
pure N-body simulations (black solid line). In the low feedback 
case, gas cools rapidly to the central regions at high z, and the 
dark matter adiabatically contracts. At z = 1 (lower panel) the 
energy transfer from gas to dark matter in the HFR has already 
considerably flattened the density profile of this latter compo- 
nent, that now clearly deviates from A^-body based expecta- 
tions. The profile of the MUGS run (LFR) is still contracted 
and has reached a logarithmic slope of a = 2. 

The creation of a core in the dark matter distribution has 
previously been attributed to rapid variations on the potential 
due to the bulk motion of gas clouds (Mashchenko et al. 2008, 
Pontzen & Governato 20 11). In Figure |4] we quantify this vari- 
ation by plotting the distance. A, between the position of the 
most bound dark matter {xom) and gas (xgas) particles. In the 
HFR (red line) the potential is rapidly changing potential is re- 
flected in the oscillations of , A with time, with the amplitude 
of the oscillations of the order of the size of the dark matter core 
(« kpc). In the LFR, the roughly constant and small value of A 
indicates a more stable potential. This indicates that the chang- 
ing potential is responsible for generating dark matter cores in 
our HFR. We note that while bulk gas motions are a natural re- 
sult of star formation and feedback, it is harder to conceive how 
such a mechanism would work with AGN feedback. 



3.1. When and how is the density profile flattened? 
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Fig. 3. — Redshift evolution of the dark matter density profile. In both 
panels the NFW fit has been obtained by fixing the concentration c,,;, = 4. 1 
and c,.,> = 6.0 at z = 4.8 and 1.0 respectively, according to results of Muiioz- 
Cuartas et al. 201 1 based on N-body cosmological simulations. 
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Fig. 4. — Evolution of the distance between the position of the dark matter 
and gas potential minima. The solid (red) line and the dashed (blue) line rep- 
resent the high and low feedback case respectively. The thin black line shows 
the Star Formation History for the high feedback run. 

Several mechanisms have been proposed to flatten dark mat- 
ter profiles, in order to reconcile the tensions between the ob- 
served cored profiles and the cusps of CDM predictions. Fully 
cosmological simulations have been able to show that such pro- 
cesses can occur in low mass (dwarf galaxy) systems (Gov- 
ernato 2010). Yet cosmological simulations of more massive 
disc galaxies have invariably resulted in adiabatic contraction 
(Gnedin et al. 2011), contradicting observed disc galaxies (Oh 
et al. 2008, Donato et al. 2009). The centrally peaked rotation 
curves, and high stellar mass fractions of these simulated galax- 



ies also fail to reproduce observations. 

In this work we have, for the first time, explicitly considered 
the feedback from radiation pressure due to massive stars in a 
cosmological hydro-dynamical simulation of galaxy formation. 
We compared the dark matter density profile of this new simu- 
lation with a twin run that only considered relatively weak feed- 
back from Supernovae. We have explicitly shown that stronger 
stellar feedback can reverse the effect of adiabatic contraction, 
and expand dark matter halos massive enough to host L* galax- 
ies . The cored profiles in our simulated galaxy have a core ra- 
dius and central density that agrees with observations (see fig. 
|2|l. Simultaneously, the high feedback simulations have rota- 
tion curves and stellar masses that are also better matches to 
observed disc galaxies than their low feedback counterparts. 

The flattening of the profile is due to the fluctuation of the 
global potential, both in its depth (Pontzen & Governato 201 1) 
and its position (fig. |4]l. The flattened profile arises at inter- 
mediate redshifts, when strong star formation and subsequent 
energy injection from feedback in shallower potential wells has 
the strongest effect. At high redshift {z ~ 3) the dark matter 
density profile is still in agreement with NFW-like predictions 
(Munoz-Cuartas et al. 2011). 

This pilot study shows that, with reasonable baryonic feed- 
back, there is enough energy input in the central region of the 
galaxy to induce rapid change in the potential and, eventually, 
induce a dark matter halo expansion, helping in reconciling ob- 
servations with CDM predictions. We emphasise that, while 
the simulation presented here does a good job in reproducing 
the observed dark matter profiles, it does not addresses all as- 
pects of the galaxy formation. 

We are working on expanding our current work to explore 
the effects of our feedback on a range of properties of galaxies 
in a wide range of masses (e.g. the two objects already shown 
in fig. |2|i. The case of Low Surface Brightness galaxies, in par- 
ticular, will be interesting since they also present cored density 
profiles yet have very low baryon fractions. Further, higher res- 
olution simulations will be required to confirm the impact of 
stellar feedback on the formation of density profiles in galaxies 
of varying masses. However, we hope this work represents a 
first step down a new path to creating more realistic galaxies. 
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